load timeseries;
N=length(timehist);
L=240; grouping=[1 5];

timehistint=linspace(0,timehist(end),N); dt=timehist(2)-timehist(1);

surfacehist=interp1(timehist,surfacehist,timehistint); surfacehist(1)=surfacehist(2);

for i=1:length(depthrange)
    qhist(:,i)=interp1(timehist,qhist(:,i),timehistint); qhist(1,i)=qhist(2,i);
    temphist(:,i)=interp1(timehist,temphist(:,i),timehistint); temphist(1,i)=temphist(2,i);
    avhist=interp1(timehist,avhist,timehistint); avhist(1)=avhist(2);
    %plot(timehist,qhist(:,i)); pause(1);
end
timehistint=timehistint/(1e6*60*60*24*365);

thermaldepth=7;
specthermal=ssa(temphist(:,thermaldepth)',L,grouping,0);
spectsurface=ssa(surfacehist,L,grouping,1);

surfacegroup=2; thermalgroup=2;
figure(3); clf; ax=plotyy(timehistint,surfacehist,timehistint,temphist(:,thermaldepth)'); set(ax(1),'ydir','rev');
figure(4); clf; ax=plotyy(timehistint,spectsurface(surfacegroup,:),timehistint,specthermal(thermalgroup,:)); set(ax(1),'ydir','rev');
figure(5); clf; subplot(2,1,1); plot(timehistint,spectsurface(surfacegroup,:)); set(gca,'ydir','rev'); subplot(2,1,2); plot(timehistint,specthermal(thermalgroup,:),'g');
figure(6); clf; hold on; plot(surfacehist-mean(surfacehist),temphist(:,thermaldepth)'-mean(temphist(:,thermaldepth)),'b.'); plot(spectsurface(surfacegroup,:),specthermal(thermalgroup,:),'r.');
figure(7); clf; ax=plotyy(timehistint,spectsurface(surfacegroup,:),timehistint,temphist(:,thermaldepth)'); set(ax(1),'ydir','rev');